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Abstract 


A three-dimensional mixed-domain PEM fuel cell model with fully-coupled transport phenomena has been developed in this paper. In this model, 
after fully justified simplifications, only one set of interfacial boundary conditions is required to connect the water content equation inside the 
membrane and the equation of the water mass fraction in the other regions. All the other conservation equations are still solved in the single-domain 
framework. Numerical results indicate that although the fully-coupled transport phenomena produce only minor effects on the overall PEM fuel 
cell performance, i.e. average current density, they impose significant effects on current distribution, net water transfer coefficient, velocity and 
density variations, and species distributions. Intricate interactions of the mass transfer across the membrane, electrochemical kinetics, density and 
velocity variations, and species distributions dictate the detailed cell performances. Therefore, for accurate PEM fuel cell modeling and simulation, 
the effects of the fully-coupled transport phenomena could not be neglected. 


© 2006 Elsevier B.V. All rights reserved. 
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1. Introduction 


Significant progress has been achieved in PEM fuel cell mod- 
eling in recent years. The existing multi-dimensional numerical 
models can generally be classified into two categories, namely 
the single-domain [1—6] and the multi-domain methods [7-10]. 
In the single-domain approach, one set of conservation equa- 
tions is established for the entire computational domain while 
special numerical treatments are applied to turn on/off a specific 
equation in specific regions of a PEM fuel cell. This renders a 
relatively simple method in terms of the mathematical formula- 
tion and numerical implementation into general CFD packages. 
However, an inherent weakness of the single-domain method lies 
in its inability to directly handle species transport, especially 
water transport, across the membrane. In order to overcome 
this deficiency and simultaneously maintain the simplicity, a 
mixed-domain approach has recently been developed [11]. In 
this mixed-domain PEM fuel cell model, the conservation equa- 
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tion of the water concentration is solved only in the gas channels, 
gas diffusion layers (GDL), and catalyst layers on the anode and 
cathode sides while a conservation equation of the water con- 
tent is solved in the membrane. These two equations are then 
connected using a set of internal boundary conditions at the 
interfaces between the membrane and the two catalyst layers. 
In addition, a fictitious water concentration treatment [12,13] 
is used for considering water transport in the membrane phase 
in the two catalyst layers. All the other conservation equations, 
including the mass, momentum, species concentrations, proton 
and electron transport, and energy, remain within the single- 
domain framework. Numerical studies have shown that this 
mixed-domain method provides accurate results in terms of not 
only the cell performances and current distributions but also the 
water content variations inside the membrane [11]. 

The mixed-domain model [11] is, however, still incomplete 
with regards to the interactions of the velocity field and the 
species transport processes. In other words, the velocity field is 
still decoupled from the species transport equations in the model. 
In order to fully integrate the various transport phenomena, the 
effects of species transport and electrochemical reactions on the 
mass and momentum conservation equations have to be taken 
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Nomenclature 


a 


See PS ol 


ETE 


xass apes 


water activity or stoichiometry parameter for 
species component 

surface area (m?) 

constant-pressure heat capacity (Jkg~! K7!) 
mass diffusivity (m? s7!) 

water content diffusivity (mol m`! s7!) 

Faraday constant, 96487C mol7! 

current density vector (A m~?) 

average current density (A m7?) 

transfer current density (A m73) 

thermal conductivity (W m -1 K7!) 

mass flux (kg m~? s7!) 

chemical components in electrochemical reac- 
tions 

stoichiometry parameter for electron 
electro-osmotic drag coefficient 

total number of species components 

net water flux through the membrane 
(molm~? s7!) 

pressure (Pa) 

interfacial flux (kg m~? s7!) 

universal gas constant (J mol`! K7!) 

source term 

temperature (K) 

velocity (m s7!) 

volume (m?) 

molecular weight (kg mol~!) 

mass fraction 


Greek letters 


a net water transfer coefficient 
ô thickness (m) 

E porosity 

@ phase potential (V) 

K proton conductivity (S m7!) 
À water content 

u chemical potential 

p density (kg m7!) 

o electronic conductivity (S m7!) 
T viscous stress tensor 
Subscripts 

cl catalyst layer 

e electrolyte or energy 

i species or interface 

m membrane, mass, or mixture 
s electron 

W water 

Superscripts 

cl catalyst layer 

eff effective value 

m membrane 


into account properly, meaning that additional source/sink terms 
have to be included in these equations, especially in the mass 
conservation equation, as proposed in [2,14,15]. This is a unique 
feature of the PEM fuel cells, distinct from the other chemically 
reacting flows, i.e. those in combustion and propulsion systems. 
This can be attributed to the fact that in PEM fuel cells, where the 
reactants are consumed is not where the products are produced. 
In fact, after the HOR reaction on the anode side the protons and 
electrons could leave the flow-field, get into the other regions, 
i.e. the membrane and the carbon material, move to the cathode 
side, react and get into the flow-field on the cathode side again. 
The processes of species moving in and out of the flow-fields, 
along with the phenomena of species being dissolved into and 
released from the membrane phase, would affect the flow-fields 
on both the anode and cathode sides. 

These phenomena have been considered in many existing 
PEM fuel cell models in the open literature [2,3,5,6,10,14]. 
However, deficiencies remain in these studies. First, models 
in the references [2,3,5,6,14] were based on the single-domain 
approach, and numerical approximations were made in handling 
water transport through the MEA. In the work of Dutta et al. [2], 
the MEA region was completely neglected from the computa- 
tional domain while in the work of Siegel et al. [3], a water 
dissolution rate in the form of the convective mass transfer was 
used, but no formulation was introduced to calculate the dissolu- 
tion coefficient. In the work of Wang and Wang [14], a fictitious 
water concentration treatment [12] was applied in MEA. As 
analyzed in [11,13], the fictitious water concentration treatment 
is theoretically valid in the catalyst layers but it lacks physical 
basis inside the membrane. It is not clear how water transport 
through the MEA was handled in [5,6]. The studies of Sivertsen 
and Djilali [10] were based on a multi-domain method with a 
fully-humidified membrane. Second, in the work of Hu et al. [5] 
and Wang and Wang [14], the conservation equations of species 
molar concentrations, instead of mass or molar fractions, have 
been solved. In reference [14], care was taken to ensure the 
consistency of the mass and species conservation equations. 

Given all these weaknesses and uncertainties in the exist- 
ing studies, the mixed-domain PEM fuel cell model [11] has 
been further extended to fully integrate the various transport 
phenomena. The formulation is based on fundamental transport 
and thermodynamic theories with assumptions properly justi- 
fied and clearly stated. This three-dimensional mixed-domain 
PEM fuel cell model is then applied to investigate the effects 
of the fully-coupled transport phenomena on the cell perfor- 
mance, current distribution, and net water transfer across the 
membrane. 


2. Theoretical formulation 


The conservation equations of mass, momentum, species 
concentrations, proton and electron transport, and energy are 
derived in the following forms 

Mass: 


V(piu) = Sm (1) 
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Momentum: 

SV (pit) = —Vpt+Vt+ Su (2) 
Species: 

V(puY;) = VoD" VY) + S; (3) 


Proton transport: 


VK" Vide) + Se = 0 (4) 
Electron transport: 

VoVo) + Ss = 0 (5) 
Energy: 

V(pCpuT) = Vk" VT) + St (6) 


An ideal equation of state is further applied herein to correlate 
pressure, density, and temperature, 


p=pRT (7) 


where the parameter, R, is the gas constant, which is related to 
the universal gas constant in the following form: 


pe (8) 
= w 
The molecular weight of a mixture, Wm, can be calculated 
as: 


N aN! 
-E : 


In the mass conservation equation, the source term, Sm, arises 
from the coupling of the flow-field and the species transport 
processes, and its exact form will be developed later in this sec- 
tion. The source terms in the momentum equations will also 
be discussed later in this section. In the species conservation 
equations, the mass fractions, Y;, are solved as the dependent 
variables. Since the summation of all the species equations will 
lead directly to the mass conservation equation, we will only 
need to solve N — | species equations, where the parameter, N, 
is the total number of the species. The source terms in these 
species equations in the catalyst layers are derived as follows: 


(a) for species except water: 


aijWi 
= (10a) 
(b) for water: 
4 iW, 
Sy = -V (ZR) Wa Iw (10b) 
nF 


In Eq. (10b), the electro-osmotic drag coefficient can be 
expressed as: 
à <14 


1.0 for aD 
Ng = 
“~ ) 480. = 14)+1.0 otherwise 


The electro-osmotic drag term at the interfaces of the mem- 
brane and the two catalyst layers will be determined using the 
interfacial boundary conditions, as developed later in this sec- 
tion. In deriving these equations, the electrochemical reactions 
are formulated in the following general form: 


X aiM; =ne” a2) 
i 


The transfer current density, j in Eqs. (10), can be calcu- 
lated using the Butler—-Volmer equation. Details are presented in 
references [16,17]. 

In this mixed-domain approach, the equation of the water 
mass fraction will be solved only in the gas channels, gas dif- 
fusion layers, and the catalyst layers on both the anode and the 
cathode sides. Inside the membrane, a conservation equation of 
the water content is established [11], 


V(D\ VA) + S = 0 (13) 


where the source term arising from the electro-osmotic drag can 
be expressed as: 


nai 
(5) (14) 


In Eq. (13), the convective effect inside the membrane has 
been neglected, and the water content diffusivity, D,, is given 
in [11]. 

A set of internal boundary conditions at the interfaces 
between the membrane and the catalyst layers on the anode 
and cathode sides is required to connect the water mass frac- 
tion and the water content equations. Based on thermodynamic 
phase equilibrium [18] and flux equality, the general forms of 
the interfacial boundary conditions are: 


S=- 


p= p” (15a) 
T™ = T” (15b) 
uP = uF (15c) 
and, 

m™ = m! (16a) 
de = 4 (16b) 
qr = 4G (16c) 


The schematic of the interface and the related boundary con- 
ditions are shown in Fig. 1. 

In Eqs. (15) and (16), the parameter i, is from 1 to N— 1. 
Furthermore, we have: 


N 
n= ai (17) 
i=1 


As described in [11], Eq. (15a) is required for solving water 
transport in the membrane once liquid water is involved, as 
suggested in Weber and Newman [19]. In the present pseudo 
single-phase PEM fuel cell model, this boundary condition is 
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Interface 


Fig. 1. Schematic of an interface and the related interfacial boundary conditions. 


neglected. Since the single-domain method is used for solv- 
ing the energy equation, Eqs. (15b) and (16b) will be closely 
approached but are not explicitly required. Eqs. (15c) and (16c) 
indicate that these interfacial boundary conditions are valid for 
all the species conservation equations. Since the concentrations 
of the other species, except water, inside the membrane are 
generally not required in PEM fuel cell modeling, the inter- 
facial boundary conditions for these species equations can be 
neglected. The cross-over fluxes of these species through the 
membrane, if needed, can be estimated approximately. There- 
fore, these species conservation equations can still be solved in 
the single-domain framework without resort to the interfacial 
boundary conditions. These cross-over fluxes are neglected in 
the present study. After these justified simplifications, only two 
boundary conditions from Eqs. (15c) and (16c) for water are 
needed at the interfaces of the membrane and the catalyst layers 
for connecting the water content equation inside the membrane 
and the water mass fraction equation in the other regions. As 
presented in [11], these two boundary conditions can be further 
developed as: 


O<a<l 


4 — J 0043 + 17.814 — 39.854? + 36.0a° 
l<a<3 


= 14+ 1.4(a — 1) 
(18a) 


z nai 
dw = Ww (-rrv. + ri) = (—pD "vy, + mYwlea 


m 


(18b) 


In Eq. (18b), water transport in the membrane phase inside 
the two catalyst layers has been considered using the fictitious 
water concentration treatment [11-13], and the effective water 
diffusivity in the catalyst layers is provided in [11]. 

After neglecting the interfacial fluxes of the other species 
except water, the mass flux across the interface in Eq. (17) 
becomes: 


m = Gy (19) 


This is the interfacial boundary condition for the mass con- 
servation equation. Since in the present PEM fuel cell model, 
velocity inside the membrane is neglected as it is very small, this 
boundary condition is required only in the two catalyst layers. In 
order to maintain the single-domain approach for the mass con- 
servation equation, this interfacial mass flux is converted into a 


source/sink term as: 
AA; 


AV, (20) 


Sm i = M x 


where the parameters, AA; and AV;, are the interfacial area and 
volume of a cell in the catalyst layer neighboring the membrane. 

Now the source/sink term in the mass conservation equation 
in the present pseudo single-phase PEM fuel cell model can be 
derived as: 


N 


Sm = X Si + Smi (21) 
i=l 


Similar expressions are also provided in [14,15]. 

In fact, the mass flux and the electro-osmotic drag term would 
impose effects on the momentum equations as well, but the 
effects are extremely small, compared with the effect of the 
Darcy’s term, and are thus neglected in the present PEM fuel 
cell model. 

The exact forms of the source/sink terms in the momentum, 
proton and electron transport, and energy equations, and the rel- 
evant physicochemical and transport parameters could be found 
in references [16,17,20]. The anisotropic transport phenomena 
in PEM fuel cells could be handled using a simplified treatment 
introduced in [21]. 

The conservation equations in Eqs. (1)—-(7) and (13), and the 
interfacial boundary conditions in Eqs. (18a) and (18b) con- 
stitute the present three-dimensional mixed-domain PEM fuel 
cell model with fully-coupled transport phenomena. The rel- 
evant source terms, physicochemical relationships, boundary 
conditions, and numerical treatments can either be found in this 
paper or follow the same procedures as in the single-domain 
approach, as detailed in references [16,17,20]. In this model, 
after fully justified simplifications, only one set of interfacial 
boundary conditions is required to connect the water content 
equation inside the membrane and the water mass fraction equa- 
tion in the other regions. All the other conservation equations 
are still solved in the single-domain framework. This model is 
applicable to the nine regions in a PEM fuel cell, including the 
current-collecting plates, gas channels, gas diffusion layers, and 
catalyst layers on both the anode and cathode sides, and the 
membrane. 


3. Result and discussion 


This numerical model has been implemented into a gen- 
eral CFD package, Fluent, and it is then applied to investigate 
the effects of the coupled transport phenomena on the detailed 
PEM fuel cell performances. The PEM fuel cell studied is a 
single straight-channel cell with a thin membrane of 25 um. Its 
geometry is illustrated in Fig. 2, and the related geometric and 
input parameters are presented in Tables | and 2, respectively. 
In the present configuration, the x-coordinate is in the through- 
membrane direction, y-coordinate the along-channel direction, 
z-coordinate the lateral direction. The cell is operated at a con- 
stant temperature of 80°C, and its inlet stoichiometry numbers 
are two based on a reference current density of 1 A cm~? on both 
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Table 1 
Cell geometric parameters 


Fuel cell geometry (mm) 


Cell length 100 
Gas channel Depth 1 
Width 1 
Layer thickness Diffusion 0.3 
Catalyst 0.01 
Membrane 0.025 
Land width 0.5 
Computational cell numbers ~320,000 


Table 2 
Input parameters 


Anode volumetric exchange current density, ajọ (A m7’) 1.0£+9 
Cathode volumetric exchange current density, ajo (A m~) 1.0£+4 
Reference hydrogen concentration, Cy, (mol m~?) 40 
Reference oxygen concentration, Co, (mol m7?) 40 

Anode transfer coefficients a= Qc =l 
Cathode transfer coefficient æc=l1 
Faraday constant, F (C mol`!) 96487 
GDL porosity 0.6 
Porosity of catalyst layer 0.12 
Volume fraction of ionomer in catalyst layer 0.4 

GDL permeability (m?) 1.0E— 12 
Equivalent weight of ionomer (kg mol™!) 1.1 

Dry membrane density (kg m~?) 1980 
Electronic conductivity in current collector (S m~!) 20000 
Effective electronic conductivity in GDL (S m7!) 5000 
Operation temperature (°C) 80 
Operation pressure (atm) 2 

Anode stoichiometry 2 
Cathode stoichiometry 2 


the anode and cathode sides. Two different cases are designed 
for the present numerical study, as shown in Table 3. Since the 
cell operates at a constant temperature of 80 °C, both cases are in 
low-humidity operation conditions, consistent with the present 
pseudo single-phase calculations. In the present numerical cal- 
culations with fully-coupled transport phenomena, hydrogen 
and water are fed into the anode inlet while air and water into the 
cathode inlet. Therefore, only the conservation equation of the 
water mass fraction is solved on the anode side, while the conser- 
vation equations of water and oxygen mass fractions are solved 
on the cathode side. Numerical results from the prior mixed- 
domain model [11], which decouples the flow-field from the 
species transport processes, are also provided for comparisons. 
The computational grid used in the present numerical study is 
the same as in the previous work [16,22], which is generated 
based on careful grid-independence studies. 


Table 3 


Inlet humidification temperature and relative humidity at cell temperature of 80°C 


MZ>nmwems 


Fig. 2. Cell geometry. 
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Fig. 3. Average current density variation along the channel direction in case 1. 


Current density variations along the flow direction from case 
1 are displayed in Fig. 3, in which the result from the present 
PEM fuel cell model with fully-coupled transport phenomena 
is compared with that from the decoupled method introduced in 
[11]. The two curves show the same trend, and the difference 
of the average current densities between them is less than 4% 
(1.13 A cm~? for the fully-coupled model versus 1.09 A cm~? 
for the decoupled model), a result agreeing with that in [14]. 
However, the maximum currents from the two models occur at 


Case number Anode 


Cathode 


Humidification temperature (°C) 


Relative humidity (%) 


Humidification temperature (°C) Relative humidity (%) 


Case 1 80 100 
Case 2 50 26 


20 5 
50 26 
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Fig. 4. Current density distribution in the mid-thickness plane of the cell inside 
the membrane in case 1: (a) from the fully-coupled model; and (b) from the 
decoupled model (unit: A m7). 


different locations; for example, the maximum current density in 
the fully-coupled model occurs in the middle of the cell length at 
y=5cm, but it shifts towards the exit in the decoupled model at 
y=8cm. The detailed current distributions from the two meth- 
ods in the mid-thickness y—z plane inside the membrane are 
illustrated in Fig. 4, in which the difference of the results from 
the two models can be clearly seen. 

Fig. 5 shows the water content distributions from the two 
models in the mid-thickness plane of the cell inside the mem- 
brane. The water content calculated from the fully-coupled 
model shows a drastic increase near the middle of the cell length, 
which contributes to the early increase of the current density, as 
shown in Fig. 4a, while the water content in the decoupled model 
varies smoothly in the along-channel direction. 

Fig. 6 compares the variations of the net water transfer coef- 
ficient in the along-channel direction calculated from the two 
models. The net water transfer coefficient is defined as the ratio 
of the net water flux through the membrane to the proton flux, 


FNw 


(22) 


a= 
Tayg 


At the first half of the cell length, the net water transfer 
coefficients from the two models are in good agreement, but sig- 
nificant difference appears starting from the middle of the cell 
length. The net water transfer coefficient from the fully-coupled 
model drops into negative values around the middle of the cell 
length but rapidly rises into positive values. It remains negative 
in the decoupled model. The difference of the results between 
the two models reaches 0.18 at the exit of the cell. The average 
value of the net water transfer coefficient from the fully-coupled 
model is positive at around 0.18, consistent with the experimen- 


13.4 

12.8 12.3 
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10.2 9.7 
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Fig. 5. Water content distribution in the mid-thickness plane of the cell inside 
the membrane in case 1: (a) from the fully-coupled model; and (b) from the 
decoupled model. 


tal analysis in [23] with a thin membrane and similar operation 
conditions. 

A simple analysis indicates that the magnitude of the trans- 
verse velocity in the through-membrane direction is proportional 
to the net water transfer coefficient in the following form: 


v ~ (Wy, + 2Wwa)l (23) 


Once the absolute value of the parameter, œ, is slightly higher 
than 0.05, it would determine the direction of the transverse 
velocity, i.e. a positive value corresponding to positive velocity 
pointing from the anode to the cathode side and vice versa. 


PO de® — fully-coupled 
de-coupled 


Along-Channel Length (cm) 


Fig. 6. Variation of the net water transfer coefficient along the channel direction 
in case 1. 
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Fig. 7. Velocity distribution in the cross section perpendicular to the flow direction in case 1: (a) at the inlet region; and (b) in the middle of the cell. 


Fig. 7 shows the velocity distributions in the x-z cross-section 
planes perpendicular to the inflow direction at two different loca- 
tions. Fig. 7a shows velocity distributions in a cross-plane near 
the cell inlet region, where from Fig. 6, the net water trans- 
fer coefficient is positive, and therefore, the overall transverse 
velocity should point from the anode to the cathode side. Fur- 
thermore, the magnitude of the transverse velocity is larger under 
the gas channel than under the land, showing the dominant effect 
of the electro-osmosis. Fig. 7b shows velocity distributions in a 
cross-section in the middle of the cell length, where according 
to Fig. 6, the net water transfer coefficient is a negative value 
at around —0.12. Consistent with the above simple theoretical 
analysis, the overall transverse velocity is now negative pointing 
from the cathode to the anode side. In addition, the magnitude of 
the velocity is larger under the land than under the gas channel, 
indicating the dominant effect of water back diffusion. 

In reference [14], a Peclet number has been used to indicate 
the relative strength of convection resulting from the transverse 
velocity to that of mass diffusion in GDL. The Pectet number is 
defined as: 


udGDL 
D 


Pe= 


(24) 


Based on the velocity magnitude displayed in Fig. 7, the 
Peclet number at the two locations in Fig. 7 is at the order of 0.03, 
indicating that the convective term is far weaker than the diffu- 
sive one. However, although this Peclet number is very small, the 
fully-coupled transport processes still show significant effects 
on the current distribution and net water transfer coefficient, as 
shown in Figs. 3, 4 and 6. Therefore, it could be concluded that 
this Peclet number in GDL is not a crucial parameter to gauge the 
significance of the fully-coupled transport phenomena on PEM 
fuel cell operations. In fact, as will be discussed in the next sec- 
tions, it is the intricate interactions of the mass transfer across the 
membrane, density and velocity variations, and species concen- 
trations on both sides of the cell that produce the drastic effects 
on current distribution and net water transfer coefficient in the 
fully-coupled PEM fuel cell model. 


Fig. 8 presents the average velocity and density variations 
along the flow direction in the anode and cathode gas channels 
scaled by their inlet values, respectively. On the anode side, 
instead of the assumed constant density and ensuing constant 
velocity along the gas channel from the flow-field calculation 
in the decoupled model, the density initially decreases, then 
increases, and finally stabilizes towards its inlet value. The veloc- 
ity steadily decreases around 60% from its inlet value. Overall, 
the anode side loses around 60% of its inflow mass through 
hydrogen consumption and water transfer across the membrane. 
In the cathode gas channel, instead of the constant density and 
velocity obtained in the decoupled model, the density gradu- 
ally decreases around 10% from the inlet to the outlet region 
while velocity increases more than 20%. On the cathode side, 
although mass is gained, density still decreases slightly. The 
reason is that it loses oxygen with higher density while it gains 
water with comparatively lower density. Overall, variations of 
the density and velocity in the cathode gas channel are small 
owing to the existence of a large amount of the inert gaseous 


Scaled Velocity / Density 


anode density 
cathode density 
anode velocity 
cathode velocity 


“0 2 4 6 8 10 
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Fig. 8. Scaled average density and velocity variations in the gas channels on 
both sides of the cell in case 1. 
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Fig. 9. Velocity distribution in the gas channel in case 1: (a) on the anode side; 
and (b) on the cathode side (unit: ms!). 


nitrogen. Detailed velocity distributions in the gas channels on 
both sides are further displayed in Fig. 9. 

Variations of the average molar concentrations of hydrogen 
and water in the anode gas channel and oxygen and water in 
the cathode gas channel along the flow direction are presented 
in Fig. 10, in which the hydrogen and oxygen concentrations 
are scaled by their respective inlet values while the water con- 
centrations on both sides are scaled by their respective outlet 
values. The molar concentration of hydrogen on the anode side 
initially increases slightly and then stabilizes towards its inlet 
value. Overall, the hydrogen concentration calculated in the 


anode H2 
=. - cathode O2 
—A— anode water 
—--A---- cathode water 


Scaled Molar Concentration 


0 2 4 6 8 10 
Along-Channel Length (cm) 


Fig. 10. Scaled average molar concentration variations in the gas channels on 
both sides of the cell in case 1. 


fully-coupled model maintains high values along the gas chan- 
nel, consistent with the experimental results in [23]. This is in 
sharp contrast with the result in the decoupled model, which 
shows steady decrease of the hydrogen concentration with the 
electrochemical consumption and is thus inconsistent with the 
experimental results [23]. The water concentration in the anode 
gas channel varies in the opposite direction of the hydrogen 
concentration; it initially decreases resulting from water loss 
dictated by the electro-osmotic drag and then increases caused 
partly by the water back diffusion and partly by the decreasing 
velocity. Towards the anode exit, water concentration increases 
slightly, although net water is lost to the cathode side in this 
region, as shown in Fig. 6. This phenomenon and the high 
hydrogen concentration values throughout the cell can both be 
attributed to the significant velocity decrease. Therefore, in the 
fully-coupled PEM fuel cell model, there exist intricate interac- 
tions between the mass transfer across the membrane, density 
and species concentration variations, and velocity changes, 
which significantly affect the detailed fuel cell performances. 

In Fig. 10, the oxygen concentration decreases almost linearly 
from the inlet to the outlet while water concentration increases 
significantly. Because the velocity and density variations on the 
cathode side are relatively small, as shown in Fig. 8, these species 
concentration variations could be easily correlated, as has been 
demonstrated in [23]. 

Fig. 11 compares the current density variations calculated 
from the fully-coupled and decoupled models in case 2. The 
two curves show the same increasing trend, and the difference 
between the average current densities is less than 4%. However, 
the maximum local difference at the exit region is very large, 
reaching 14%. The variations of the net water transfer coefficient 
from the two models are shown in Fig. 12. The values obtained 
from the present PEM fuel cell model with fully-coupled trans- 
port phenomena are always higher than those obtained from the 
decoupled model, and the gap becomes larger towards the exit 
region. At the exit, the gap reaches 0.06, showing more than 
26% relative difference. 

Variations of the scaled average hydrogen and water con- 
centrations in the anode gas channel and oxygen and water 
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Fig. 11. Average current density variation along the channel direction in case 2. 
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Fig. 12. Variation of the net water transfer coefficient along the channel direction 
in case 2. 
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Fig. 13. Scaled average molar concentration variations in the gas channels on 
both sides of the cell in case 2. 


concentrations in the cathode gas channel in case 2 are pre- 
sented in Fig. 13. The hydrogen concentration on the anode side 
decreases slightly, but as in case 1, it maintains very high value 
throughout the entire cell. 


4. Conclusion 


In this paper, a previously developed three-dimensional 
mixed-domain PEM fuel cell model has been further extended 
to fully integrate the various transport phenomena. In this 
model, after fully justified simplifications, only one set of 
interfacial boundary conditions is required to connect the water 
content equation inside the membrane and the water mass 


fraction equation in the other regions. All the other conservation 
equations are still solved in the single-domain framework. This 
numerical model is applicable to the nine regions in a PEM 
fuel cell, including the current-collecting plates, gas channels, 
gas diffusion layers, and catalyst layers on both the anode and 
cathode sides, and the membrane. 

Although the fully-coupled transport phenomena produce 
only minor effects on the overall PEM fuel cell performance, i.e. 
average current density, they make significant effects on current 
distribution, net water transfer coefficient, velocity and density 
variations, and species distributions in a PEM fuel cell. There- 
fore, for accurate PEM fuel cell modeling, the fully-coupled 
transport phenomena could not be neglected. 

Numerical results indicate that the Peclet number defined 
in GDL, measuring the relative strength of convection to mass 
diffusion processes in the region, is not a crucial parameter to 
gauge the effects of the fully-coupled transport phenomena on 
the PEM fuel cell operations. In fact, there exist intricate interac- 
tions of the mass transfer across the membrane, electrochemical 
kinetics, density and velocity variations, and species distribu- 
tions, which play a dictating role in determining the detailed 
cell performances. 
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